function expandFromPoint2(inputFileName,sliceNumber, seedDiameter, outputFileName)
close all;
%% Load and display the file

origImage = load_nifti(inputFileName);

targetImage = origImage(:,:,sliceNumber);
himage = Imshow(targetImage);
[startXX,startY] = ginput(1);
% convert the point to int
startX = int8(startY)
startY = int8(startXX)

%%
% define additional params
thresVal = double((max(targetImage(:)) - min(targetImage(:)))) * 0.05;
[rowNum, colNum, slidesNum] = size(origImage);

%TODO remove
slidesNum = 3;
targetPixelVal = double(origImage(startX, startY, sliceNumber));

% Init the result matrix (Seeded Region Growing)
SRG = false(size(origImage)); %false(n) is much faster and more memory efficient than logical(zeros(n)).

%% Fill the queue with start points (Start point +/- seedDiameter)
queue = [];

for i = -seedDiameter:seedDiameter
    for j = -seedDiameter:seedDiameter
        for k = sliceNumber-seedDiameter: sliceNumber +seedDiameter
            isPointInside=(startX +i >=1)&&(startY +j >=1)&&(sliceNumber +k >=1)&&(startX +i<=rowNum)&&(startY +j<=colNum)&&(sliceNumber+k<=slidesNum);
            if  isPointInside
                queue(end+1,:) = [startX+i, startY+j, sliceNumber+k];
            end
        end
    end
end

%% Start

while size(queue, 1)
    % pop the first element from queue (BFS)
    xv = queue(1,1);
    yv = queue(1,2);
    zv = queue(1,3);
    
    queue(1,:) = [];

    % check the neighbors for the current position
    for i = -1:1
        for j = -1:1
            for k = -1:1
                if  any([i, j, k])      % [0,0,0] => root point (already checked)
                    isNeighbourInside=(xv +i >=1)&&(yv +j >=1)&&(zv +k >=1)&&(xv +i<=rowNum)&&(yv +j<=colNum)&&(zv+k<=slidesNum);
                    if(~isNeighbourInside)
                        continue
                    end
                    isPositionSeted = SRG(xv+i, yv+j, zv+k); %this pixel checked already
                    currentStd = std(stdMatrix)
                    origImage(xv+i, yv+j, zv+k)
                    
                    isInRange = origImage(xv+i, yv+j, zv+k) <= (targetPixelVal + thresVal) &&...
                                origImage(xv+i, yv+j, zv+k) >= (targetPixelVal-thresVal);
                            
                    if ~isPositionSeted && isInRange

                        queue(end+1,:) = [xv+i, yv+j, zv+k]; %add neighbour to BFS queue
                        SRG(xv+i, yv+j, zv+k) = true; % connect the pixel to region
                        % TODO
                        targetPixelVal = mean(mean(origImage(SRG > 0)));

                    end
                end
            end
        end
    end
end

%Finish




%%
figure();
Imshow(SRG(:,:,2));

end